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We consider solvability of the generalized reaction-diffusion equation with both space- and time- 
dependent diffusion and reaction terms by means of the similarity method. By introducing the 
similarity variable, the reaction-diffusion equation is reduced to an ordinary differential equation. 
Matching the resulting ordinary differential equation with known exactly solvable equations, one can 
obtain corresponding exactly solvable reaction-diffusion systems. Several representative examples 
of exactly solvable reaction-diffusion equations are presented. 
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I. INTRODUCTION 

Many natural phenomena involve the change of concentration/population of one or more substances/species dis¬ 
tributed in space under the influence of two processes: local reaction which modify the concentration/population, 
and diffusion which causes the substances/species to spread in space. Such phenomena are well modelled by the 
reaction-diffusion equation (RDE). 

The general form of RDE of the concentration W(x, t) of a single component in one spatial dimension is 

-W{x,t)=D—W{x,t) + f{W), (1) 

where D is the constant diffusion coefficient and f{W) is the reaction term which accounts for the local reaction. 
Eq. ((T|) is also called the KPP (Kolmogorov-Petrovsky-Piscunov) equation, named after the authors who first studied 
some of the mathematical properties of the RDE. This equation encompasses the diffusion (heat) equation (/ = 0) 
and the Fokker-Planck equation (when / is a gradient term of some function linear in W) [^. 

Different forms of the reaction term / have been proposed to describe different phenomena. For instance, the 
choice f = W{1 — W) yields the Fisher equation employed in the study of wave propagation of advantageous genes 
in a population Q and evolution of a neutron population in a nuclear reactor [^. Rayleigh-Bernard convection 
is studied using RDE with / = WQ — W^) i: while combustion and shock waves phenomena invoke RDE with 
f = W{1 — W){W — a)(0 < a < 1) Q. Generalization and extension of the KPP equation to higher dimensions and 
multi-component cases also find interesting applications in chemical kinetics Q, pattern formation and morphogenesis 
0 , nerve pulse propagation in nerve systems , and other biological systems [l^ . 

In view of its broad applicability, it is thus desirable to obtain analytic solutions of the RDE for as many systems 
as possible. However, just as any equation in sciences, solving the RDE exactly is in general a formidable task, except 
in a few simplified cases. Fortunately, for many cases of the RDE mentioned before, exact solutions can be found in 
the form of travelling wave solutions |lll| . 

In this paper we would like to consider exact solvability of the RDE in terms of the similarity solutions [l^ • This 
is motivated by our recent works on similarity solutions of the Fokker-Planck equation, which as mentioned before is 
a subclass of the RDE [l3l - [l^ . We found it very interesting that for a class of the Fokker-Planck equation with time- 
and space-dependent coefficients, a general formula of exact solutions can be obtained in closed form by the similarity 
method, both for fixed and moving boundaries. One advantage of the similarity method is that it allows one to reduce 
the partial differential equation under consideration to an ordinary differential equation which is generally easier to 
solve, provided that the original equation possesses proper scaling property under certain scaling transformation of 
the basic variables. Here we would like to extend our previous consideration to the RDE. It turns out that, similar 
to the Fokker-Planck case, one can determine the possible functional forms of the diffusion and the reaction term of 
the RDE in order to get exactly solvable system with similarity solutions. 

This paper is organized as follows. Sect. H discusses the scaling properties of the RDE. Sect. HI introduces the 
corresponding similarity variable and scaling forms of the relevant functions, which are used to reduce the RDE into 
an ordinary differential equation. The equation of continuity is discussed in Sect. IV which helps to identify two types 
of scaling behaviours of the RDE. Some examples of these two types of scaling reaction-diffusion (RD) systems are 
presented in Sect. V and VI, respectively. Sect. VH concludes the paper. 
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II. SCALING OF REACTION-DIFFUSION EQUATION 


We shall consider the following general form of the RDE in (1 + l)-dimension 

dW{x,t) 


dt 




( 2 ) 


where W{x,t) is the particle density function, D{W,x,t) is the diffusion coefficient and f{W,x,t) the reaction term. 
We use the term “particle” to denote generally the number of basic member of a substance or a specie. The domains 
we shall consider in this paper are the real line x G (— 00 , 00 ), or the half lines x G [0,oo) and (—ooffi]. Cases with 
finite domains, which correspond to systems with moving boundaries, can be considered similarly |14| . To cater for 
the most general situation, we leave the possibility that D and / could be functions of W. 

We shall consider the similarity solutions of the RDE. Such solutions are possible, provided the RDE possesses 
certain scaling symmetry. Below we shall study the scaling property of the RDE. 

Consider the scale transformation 


X = e“x , t = 


( 3 ) 


where the scale factor e and the two scaling exponents a and b are real parameters. Suppose under this transformation, 
the density function, the diffusion coefficient and the reaction term scale as 

W{x,t) = e^W{x,t), D{W,x,t) = e'^D{W ,x,t), f{W,x,t) = e‘^f{W,x,t). (4) 

Here the scaling exponents c, d and e are also some real parameters. Written in the transformed variables, Eq.([2]) 
becomes 




( 5 ) 


For simplicity and clarity of presentation, here and below we shall often omit the independent variables in a function. 

One sees that if the scaling indices satisfy b = 2a — d = c — e, then Eq.([5]) has the same functional form as Eq.(l2|). 
In this case, the RDE admits similarity solutions. We shall present such solutions below. 


III. SIMILARITY VARIABLE AND SCALING FORMS 

The similarity method is a very useful method for solving a partial differential equation which possesses proper 
scaling behavior. One advantage of the similarity method is to reduce the order of a partial differential equation 
through some new independent variables (called similarity variables), which are certain combinations of the old inde¬ 
pendent variables such that they are scaling invariant, i.e., no appearance of parameter e, as a scaling transformation 
is performed. 

In our case, the second order RDE can be transformed into an ordinary differential equation which is generally 
easier to solve. Here there is only one similarity variable z, which can be defined as 

X CL 

z = —, where a = - , b ^ 0 . (6) 

Next we assume the following scaling forms of the density function, the diffusion and the reaction terms in terms 
of z: 


W{x,t) =t^y{z), D{W,x,t) =t’'p{z), f{W,x,t) =t^a{z). 

Erom Eq. 0 together with the scaling conditions b = 2a — d = c — e, one has 

c d e 

^=-, v=- = 2a-l, \=-=p-l. 

Thus a and p are the only two independent scaling exponents of the RDE. 

In terms of these scaling forms, Eq. Q reduces to an ordinary differential equation 


d 

dz 



az- - py + cryz) = 0. 

dy 


( 7 ) 

( 8 ) 


( 9 ) 
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Note that when = —a, which we will encounter below, Eq. ([9]) reduces to 


d 

dz 


p-^y + oizy ) + a[z) = 0. 


( 10 ) 


To proceed further, we shall consider the conditions imposed by the continuity in the change of the particle number 
of the system, i.e., the equation of continuity. 


IV. EQUATION OF CONTINUITY 

The total number N of the system is related to the density function W{x,t) by 

N = [ W{x,t)dx = t°‘+f^ [ y{z)dz, 

Jv Jt> 


( 11 ) 


where T) is the domain of the independent variable. For simplicity, we use the same notation T? for both the variable 
X, and the corresponding similarity variable z. 

Eq.d]) distinguishes two different situations: a + y ^ 0 and a + p, = 0. It is obvious from this equation that N is 
conserved if and only if y = —a. 

Eq.dlll) implies that 


dN 

dt 


= (a + /r)t“+^-i (^J^yiz)dzy 


( 12 ) 


On the other hand, from Eq. ([2]) one has 
dN 
dt 


= J ^^dx = AyD{W,x,t)y-W{x,t)^ +J f{W,x,t)dx 

= [ a{z)dz. 

\ dz J Jx> 


(13) 


Here dV denotes the boundaries of the domain 2?, and A(- • • )ax> the difference of the terms in the bracket at the 
boundaries. 

In view of A = /i — 1, one has 


{a + p)J y{z)dz = j a{z) dz + A 

If y{z) —0 fast enough so that the boundary terms tend to zero, then one has 

{a + p.) / y{z)dz= / a{z)dz. 

Jv Jv 

This is the situation we shall consider in most of the cases below, except the example in Sect. VI.D. 
In what follows, we shall present some examples for p = —a and p ^ —a. 


(14) 


(15) 


V. CASES WITH p = -a 

As mentioned at the beginning of Sect. IV, N is conserved when p = —a. Hence one can normalize W (x, t) = t~°‘y{z) 
and consider it as the probability distribution function. 

Furthermore, Eq. m implies a{z)dz = 0. This is most easily satisfied if cr(z) is a total differential, i.e., 
a{z) = —dT{z)/dz for some function t{z). This includes as subclass any function cr(z) that is anti-symmetric w.r.t. 
the mid-point of the domain D in the similarity variable z. But this latter situation is possible only if D is the whole 
line or a finite domains in the z-space (corresponding to moving boundaries in the x-space), and is not possible for 
the half-line. 

Consequently, Eq. m is a total derivative, and can be integrated once to give 

py' -I- azy — r = constant. (16) 

Here the prime denotes derivative w.r.t z. For y, r —)• 0 at the boundaries, the constant equals zero and we need only 
to consider the following equation instead 


py' + azy — r = 0. 


(17) 
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A. Fokker-Planck type 


First we consider the situation in which the function t{z) is proportional to y{z), i.e. t{z) = f3[z)y{z) for some 
function /3(z). In this case the RDE is of the Fokker-Planck type, where the function /3(z) plays the role of the drift 
coefficient. Integrating Eq. (II3 once gives 


y(z)ocexp^y dz . (18) 

Hence for any choice of I3{z) and p{z) such that y{z) in Eq. (fTSl) is integrable and that W{x, t) is normalizable, one has 
an exactly solvable RD system. This is exactly the same as the way to obtain similarity solutions of the Fokker-Planck 
equations discussed in [la - ll^ . All the cases presented there for the Fokker-Planck equations can be carried over to 
this type of RDE. As such we shall be brief on this case, and present only an example for illustration. 

Let us take p{z) = I and (i{z) = /3iz + /3o {/3i < a). Then the system is given by D{x,t) = , f{x,t) = 

and 


W{x,t) = 


a-Pi -1 


-Pi 






a > Pi, 


(19) 


for t > 0 and —oo < x < oo. This solution represents a diffusing wave wave with a moving peak. 

Other possible forms of y(z), including those related to solutions with moving boundaries, and solutions involving 
the recently discovered exceptional orthogonal polynomials, can be found in Ref. 141 and [1^, respectively. 


B. Non-Fokker-Planck type 

We now consider the situation where a(z) is not proportional to y(z), i.e., cr(z) ^ y[z). 

Again we assume ?/,r —> 0 at the boundaries, so the constant in Eq. (1161) equals zero and we need consider only 
Eq. (fT71) . The general solution of Eq. (IT71) is 

(17 = constant. (20) 

Any choice of p{z) and t{z) such that y{z) is exactly integrable and W{x,t) is nomalizable furnishes a solvable RD 
system. 

Some examples below serve to elucidate the idea. 


1. p{z) = 1 

Take t{z) = 27 zexp(—ryz^) with real constants 7,77 > 0. The general solution y{z) in Ea. (l20l) is 

y{z) = -b ( 21 ) 

a — IT] 

for z G (— 00 , cx)). Thus we obtain an exactly solvable RD system with 

D{x,t) = 

f{x,t) = -27 y - 2r] ( 22 ) 

W{x,t) = t~° ( 

\a-2r] J 


To ensure W{x, t) > 0 for all x £ (— 00 , 00 ), we must have C > =F 27 /|a — 277 I for a > 2p and a < 2r], respectively. 
In Eig.l we show the graphs of D{x,t), f{x,t) and W{x,t) for a set of parameters with three different times. 
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2. p{z) = z 

Let t(z) = 2 exp(— 772 ) with real constants 77 > 0. The general solution y{z) is 

y{z) = + C7e-“^ 

a — Tf 


(23) 


for 2 G [0, 00 ). 

The corresponding exactly solvable RD system is defiined by 

D{x,t) = 

f[x,t) = (24) 

W{x,t) = r°‘ . 

\a-y J 

To ensure W{x,t) > 0 for all x G [0, 00 ), we must have C > Tl/|ci^ — ‘n\ for a > ry and a < y, respectively. 

We show in Fig. 2 the graphs of D{x^ t), f{x, t) and VF(x, t) for a set of parameters with three different times. 


3. p{z) = az^ 

For p{z) = az"^, Eq. (|20|) becomes 

*'<*> = i (/ ^ 

One must set C = 0 in order that y{z) is finite at z = 0. 

As an example, let us take t(z) such that y{z) = exp(—/3z^),/3 > 0. This means 

t{z) = az (1 — 2/3z^) . 

The corresponding RD system is defined by 

1 


(25) 


(26) 


D{x,t) = t 0 < a; < 00 , 


f[x,t) = 




- 1 




(27) 


W{x,t) = 


It is interesting to see that in a rather complicated reaction term /(x, t) leads to a simple Gaussian-like distribution 
W{x,t). 


C. A nonlinear diffusion equation 

Taking t(z) = 0 and p{z) = y"'{z) (n = 0,1, 2,...) , one gets a nonlinear equation 


Iz «2?/(z) ) = 0. 


(28) 


This is the reduced scaling ODE associated with the nonlinear diffusion equation 


where the diffusion coefficient is given by D(x,t) = kF"(x,t) [Tf 


(29) 
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In [l^ . Eq. (1291) was solved through a very indirect nonlinear transformation of the function W{x,t), namely, 

w 

f^dt. (30) 

We shall see that the same solution was extremely easy to derive based on the similarity method here. 

First we note that the form oi D = implies a relation between the scaling exponents d = nc, which leads to 
the relation v = n^j,. From v = 2a — 1, one gets 

a = -/r=^—, n = 0,1,2,... (31) 

n + 2 

Next by integrating Eq. (1^ once, one gets 

?/"(z)-^ 2 /(z) + azy{z) = constant (32) 

dz 

The case with constant 0 is rather complicated, so here we shall follow [I^ and consider only the case with 
constant = 0. In this case, Eq. (l32l) is easily integrated to give 



y[zj = - —z J , C = constant. 

From W{x,t) = t^y(z) and y = —a = —l/(n + 2), we arrive easily at the solution given in 


(33) 


W{x,t) 


cr 


2{n + 2) t 


(34) 


(35) 


VI. CASES WITH ^ 7 ^ -a 

Now we come to the second type of scaling behaviour of the RDE, with y ^ —a. For such RD system, the number 
of particles does not conserve. 

The relevant equation is Eq. m 


py" + {p' + a z)y' + a- yy = Q 


(36) 


A. p = 1 


In this case, Eq. dHI becomes 


y"{z) + azy'{z) - yy{z) + (j{y) = 0 . 


(37) 


One can match this equation with any equation whose solutions are exactly known. This will then determine the 
function <j{y) and all the parameters. 

For instance, the differential equation in §2.1.2.29 of 


y" + {az + h)y' + c[{a - c)z‘^ + hz + l]y = 0, 
is known to have a particular solution 

2/0 = e 2“ , 

for arbitrary constant c. For our purpose, we must take c > 0 to ensure normalizability of W{x^t). 
Comparing the two equations, one finds that for the choice 

a = a, 6 = 0, c > 0, cj{y) = [c{a — c)z^ + {p + c)]yo{z), z G (—oo, oo). 


(38) 


(39) 


(40) 
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we obtain an exactly solvable RD system with 

D{x,t) = 2 /(z) = ^ xe(— 00 , 00 ), 

f{x,t) = c(a-c)(^) +(m + c) , (41) 

W{x,t) = 

It can be checked that the continuity equation (fT^ is satisfied. 

In Fig.3 we show the graphs of D{x,t), f{x,t) and W{x,t) for a set of parameters with three different times. 


B. p = P z 


The relevant equation is 


zy”{z) + (^^z + 1^ y'{z) + 


iy) - yyiz) 


= 0 . 


This can be matched with the equation in §2.1.2.71 of [13 

zy" + {az + b)y' + c[(a — c)z + b\y = 0, 


(42) 


(43) 


which has a particular solution yo = e “ for z G [0, 00 ). We shall take c > 0 so that y{z) ^ 0 as a; —)■ 00 . To match 
Eq. (1421) with this equation, one can take a = a//3, 6 = 1, and 


Then 


a{y) = [c(a - /3 c)z + y. + (5 c] yo- 

D{x,t) = pt°‘~^x, y(z) = a;G[0,oo), 

f{x,t) = \c{a-Pc) ^ + {y +Pc) 


--C-^ 
g 2 ^ 


(44) 


(45) 


W{x,t) = t^e 2 ‘=t° . 


C. p=-zaz^ 


Next we consider the choice p = —^az^. This choice eliminates the first derivative term in y{z). The RDE becomes 


y — 


2 f cr - yy 


= 0 . 


a \ 

With the choice a{y) = [^az'^{g'^ + g') + p\ y{z), Eq. (l46)l reduces to 

y” - {9^ + g')y = O- 

This equation admits a particular solution 


(46) 


(47) 


(48) 


We shall choose g{z) such that y{z) —>■ 0 fast enough that the boundary terms go to zero. For the above choice of 
a{y) and yo(z), we obtain an exactly solvable RD system with 


1 

D{x,t) = 


f{x,t) = ^cr{y), 

W{x,t) = exp g{z)d^ . 


(49) 


Thus an exactly solvable RD system is obtained for every choice of g{z) such that it is integrable. Below we shall 
present two simple examples. 








1. 


3 = -l 


In this case we have 


yo{z) 


e , 




W{x,t) = t^e 


for t > 0 and x S [0, oo). 


(50) 


a. g = -z 


In this case we have 


yo{z) 

W{x,t) 


-'{Kf) [(f)' 


-1 




(51) 


for f > 0 and x € (—oo, oo). 


D. Nonlinear case: Generalized Fisher equation 


Let us consider the situation where the coefficients of Eq. (l36l) of y" and y' are proportional, say, 

p' + az = 7 p. 


The general solution of this equation is 


P{z) 


Q/ 1 

— (z -\—) + /3e^^, (3 = constant. 
7 7 


(52) 


(53) 


Eq. (1551) then reduces to 


y” + iy' + -{(j - py) =Q. (54) 

P 

An exactly solvable RDE can be obtained if one can match this equation with a known solvable ODE of the same 
form. However, it can be easily checked that, to satisfy the continuity equation (I14L the parameter a must be zero, 
i.e. a = 0, and thus z = x,v = —1 and D{x, t) = 

As an example, let us consider the generalized Fisher equation 


y" + jy' + y{l-y'^) = 0, n > 0, (55) 

which is an important equation in mathematical biology and nuclear physics ii®. To match the RDE requires 

(7 = /ry + (1 - y”). (56) 


The corresponding RDE is 


(57) 


This equation is invariant under the scale transformation 

x = x, t = eH, W = 


(58) 
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It is known that particular solutions of Eq. (l55l) are possible if 7 = ±(h„ + l/h„), where = Y^n/2 + 1 (Tsl - f^ . 
The particular solutions are 

2 /±(z) = + (59) 

where C is some real constant. We will take C > 0 so that W{x, t) is regular. Also, the system defined by 7 < 0 
is just the mirror image of that with 7 > 0, as Eq. (15^ is invariant under 7 —?► —7 and x —>■ —x. So it suffices to 
consider just the case with 7 > 0 . 

In Fig. 4 we show the graphs of D{x,t), f{x,t) and W{x,t) for a set of parameters with three different times and 
negative exponent ^ < 0. D{x,t) is independent of /i, and f{x,t) and W{x,t) decrease as time increases. For positive 
/i > 0, both f{x,t) and W{x,t) will increase as time increases. The traveling wave behaviour of the solution of the 
generalised Fisher equation is here transformed into a scaling wave behaviour. 


VII. SUMMARY 

We have considered solvability of the reaction-diffusion equation with both space- and time-dependent diffusion 
and reaction terms by means of the similarity method. By introducing the similarity variable, the reaction-diffusion 
equation is reduced to an ordinary differential equation. It is interesting to realise that the reduced ordinary differential 
equations, namely, Eqs. ([9|) and (flTl) . are quite simple in their functional forms. Particularly, Eq. (ED is integrable and 
its solution can be given in closed form. By matching these two ordinary differential equations with known exactly 
solvable equations, one can obtain corresponding exactly solvable reaction-diffusion systems. We have presented 
several representative examples of exactly solvable reaction-diffusion equation. 

Of course, similarity solutions, just as the travelling wave solutions, are only one type of many possible symmetry 
solutions one can consider for the RDF. Under different conditions, the RDF may admit other symmetry solutions. 
Very recently, classification of exactly solvable RDF with gradient-dependent diffusivity (for D a derivative of W, 
i.e., D = D{Wx)) has been considered using the Lie symmetry approach [2l|. It would be of interest to extend such 
consideration to more general diffusion and reaction terms. 
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FIG. 1: Plot of D{x, t), f(x, t) and W(x, t) for a = 0.6, 7 = 2,77 = 0.1, G =1 and time t = 1.0 (dotted), t = 2.0 (dashed), t = 3.0 
(solid). 





FIG. 2: Plot of D{x,t), f{x,t) and W{x,t) for q = 2,?7 = 1,G=1 and time t = 1.0 (dotted), t = 1.25 (dashed), t = 1.5 (solid). 



FIG. 3: Plot of D{x, t), f{x, t) and W{x, t) for a = 1, p = 0.5, C = 1 and time t = 1.0 (dotted), t = 2.0 (dashed), t = 3.0 (solid). 
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FIG. 4: Plot of D{x,t), f{x,t) and W{x,t) for n = 3,/r = —1,13 = 1, C = 1 and time t — 0.05 (dotted), 0.1 (dashed) and 0.2 
(solid). 





















